# 4.8

rm(list=ls())
library(ggplot2)
library(latex2exp)



x <- c(25669,17885,32070,13050,18128,22247,
       14777,15386,28179,77388,47251,24408,
       28811,18499,68024,40472,32665,31551,
       80855,18318,4053,17741,32935,11777,
       14788,1151,19400,7200,2572,3169,9650)
y <- c(7944,12779,31750,14198,15080,6692,13923,
       10648,6756,49663,30276,27033,23237,19694,
       53323,40415,30012,28353,33304,18237,3890,
       16048,28812,13204,16119,1596,20825,9664,
       3528,3794,10288)

model = lm(y~x)
summary(model)
cor.test(x,abs(model$residuals),method = "spearman")

# WLS

m = c(-2,-1.5,-1,-0.5,0,0.5,1,1.5,2,2.5,3,3.5,4,4.5,5)

LOGLIKE = - Inf
best_m = -2
for(i in m)
{
  model_temp = lm(y~x, weights = (x)^{-i})
  LOGLIKE_c = logLik(model_temp)
  if(LOGLIKE_c>LOGLIKE)
  {
    LOGLIKE = LOGLIKE_c
    best_m = i
  }
  print(i)
  print(LOGLIKE_c)
}


model_wls = lm(y~x, weights = (x)^{-2.5})


data_1 = data.frame(x = x, y = y)
data_2 = data.frame(x = x, y_hat = model$fitted.values,method=rep('OLS',31))
data_3 = data.frame(x = x, y_hat = model_wls$fitted.values,method=rep('WLS',31))
data_23 = rbind(data_2,data_3)

p <-  ggplot(data = data_23)+ 
      geom_line(aes(x = x, y = y_hat, group = method,color=method))+
      geom_point(data = data_1, aes(x = x, y = y))+
      labs(x="x", y= "y")+
      theme_classic()

ggsave('fit_4_8.png',p, width=6,height=5)
p

